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Abstract 

We review recent molecular dynamics simulations of thermally activated undulations and defects in the lamellar La phase of 
a binary amphiphile-solvent mixture, using an idealized molecular coarse-grained model: Solvent particles are represented by 
beads, and amphiphiles by bead-and-spring tetramers. We find that our results are in excellent agreement with the predictions 
of simple mesoscopic theories: An effective interface model for the undulations, and a line tension model for the (pore) 
defects. We calculate the binding rigidity and the compressibility modulus of the lamellar stack as well as the line tension 
of the pore rim. Finally, we discuss implications for polymer-membrane systems. 
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1. Introduction 

Amphiphilic molecules, such as lipids, contain hy- 
drophilic (water loving) and hydrophobic (water fear- 
ing) parts. In water, they self-assemble spontaneously 
such that the hydrophobic parts are shielded from the 
water environment. One of the most prominent struc- 
tures is the lamellar phase, where the amphiphiles are 
arranged in stacks of bilayers. Lipid bilayers are of 
special interest, because they are essential components 
of biological membranes. 

Phenomenologically, membranes are often de- 
scribed as two dimensional surfaces. A stack of 
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nearly planar membranes of average distance d can 
be parametrized by a set of functions m,i(x), where 
X = {xi,X2) are coordinates in the plane and x-^ — 
Un + nd denotes the position of the nth membrane in 
the direction perpendicular to the plane. The simplest 
theoretical Ansatz approximates the free energy of 
such a stack by [1] 



U„-(_i— u„ 



(1) 



with the bending rigidity and the compressibil- 
ity B. This free energy functional describes thermally 
activated membrane undulations with a characteristic 
in-plane correlation length £, — [KclB)^l^. 

Apart from undulations, membrane stacks also 
contain defects. Among these, the membrane pores 
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are particularly important, because they influence 
crucially the permeation properties of the membrane. 
In the simplest phenomenological model for pore 
formation [2], a pore with the area a and the contour 
length c has the free energy 



-Epore = -Bo + Ac - 7a, 



(2) 



where 7 is the surface tension of the membrane and 
A a line tension. In an isotropic membrane stack, the 
membranes are tensionless with 7 = 0. 

In the present paper, we describe large scale molec- 
ular dynamics simulations of a molecular coarse-grain 
model which allow to test these simple theoretical con- 
cepts. We find that the models models (1) and (2) cap- 
ture the physics of our membrane stacks very well. 
The results are then used to discuss a more complex 
system, i.e., a polymer inserted in a membrane stack. 
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Fig. 1. Total density profiles (top) and partial volume fractions 
of head, tail and solvent beads. (Head and tail beads refer to P 
and H beads of amphiphiles). 

3. Results 



Membrane Structure and Thermal Fluctuations 



2. Simulation Model, Method, and Data Analysis 

The simulation model was derived from a similar 
model by Soddemann et al [3]. The "molecules" con- 
sist of two types of beads, H or P. "Solvent" par- 
ticles are single P beads, "Amphiphiles" are H2P2 
tetramers. All beads have the same mass m and diam- 
eter a and interact with a repulsive hard core potential. 
In addition, beads of the same type attract each other. 
At appropriate parameter values, the amphiphiles self- 
assemble spontaneously and form a lamellar stack of 
bUayers. 

Our simulations were carried out in the NPT en- 
semble with a Langevin thermostat. All sides of the 
simulation box were aUowed to fluctuate in order to 
ensure isotropic pressure. We studied systems contain- 
ing a total number of 153600 beads, 30720 of which 
were solvent particles, at pressure P = 2.9kBT/a^. 
The total run length was lO'^ MD steps, corresponding 
to the time span 10^ ^rna'^ /ksT. The initial equili- 
bration time was 10^ MD steps. Further details can be 
found in Refs. [4,5]. 

To analyze the data, the simulation box was subdi- 
vided into columns of diameter ~ 1.3a, perpendicular 
to the membranes. For each configuration and in each 
column, the local profiles were determined. This al- 
lowed to locahze the membranes and to identify mem- 
brane defects [4,5]. 



The local membrane structure can be characterized 
by the local density and composition profiles. In order 
to eliminate the smearing effect of the undulations, we 
calculate the profiles in each column separately and 
shift it by the local interface position before averaging. 
The results are shown in Fig. 1 . The total density varies 
very little throughout the system, but the amphiphile 
and solvent layers are nevertheless well-separated. 

From the statistics of the interface positions, we 
can analyze the capillary wave spectrum. Here we 
show only the transmembrane structure factor be- 
tween membranes, which we define as Sn{q) = 
(■u„(q)*uo(q)), where u„(q) is the two dimensional 
Fourier transform of u„(x). Within the model (1), 
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Fig. 2. Transmembrane structure factor, fitted to the prediction of 
the model 1 with one fit parameter ^ = {Kc/By/*. 
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we can calculate the quantity Sn{q) theoretically and 
compare the results with the simulation data. The 
only fit parameter is the parallel correlation length 
^ = {Kc/B)^^/^\ Fig. 2 shows that the theoretical 
prediction fits the data very well at ^ = 2.34cr. 

Similar fits allow to determine the layer dis- 
tance d = 6.38(7, and the Caille parameter rji = 
7:kBT/2(P^K^ = 0.053. This gives us all the 
parameters of the model (1), i.e., Kc — ^ksT and 
B = 0.13A;bT/o-4. 

Pore Defects 

The main defects in our system were pores, and 
only these were studied in detail. The density pro- 
files of the pores (not shown here, see Ref. [5]) show 
unambiguously that larger pores are hydrophilic, i.e., 
the amphiphiles at the pore rim rearrange themselves 
such that the solvent in the pore is shielded from the 
hydrophobic molecules tails. Since our amphiphiles 
are very short, the energy penalty on pores is low and 
the number of pores in the membranes is compara- 
tively large. Nevertheless, we have established in sev- 
eral tests that the pores are essentially uncorrelated. 
Therefore, pore interactions can be neglected, and we 
can attempt to interpret our data in terms of the model 
(2) with 7 = 0. 

The validity of Eq. 2 was tested in several ways. 
First, we verified that the pore contour lengths c were 
distributed according to a Boltzmann distribution. This 
allowed to extract the value of the line tension, A = 
O.lksT/a. Next we studied the shape of the pores. 
Since the membranes have no surface tension, the 
pores are not circular: At 7 = 0, Eq. 2 describes an 
ensemble of self avoiding ring chains. These should 
have by the fractal dimension of two-dimensional self- 
avoiding polymers, characterized by the Flory expo- 
nent f2 = 3/4. For example, the area a of the pores 
scales with the contour length like a oc c^/^. Fig. 3 
a) shows that the simulation data are compatible with 
this scaUng behavior. 

The energy model 2 also has implications for the 
pore dynamics: The dynamical evolution of the pore 
contour lengths should correspond to a random walk in 
a linear potential, E = Eo + Ac. This stochastic model 
has been solved analytically [6] and we can compare 
our data with exact results, e.g., for first passage time 
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Fig. 3. Tests of the line tension model (2). Top: Mean area a of 
pores as a function of contour length c. Bottom: Probability that a 
pore with contour length c = 10 reaches the contour c = 8 for the 
first time at the time t, compared with the theoretical prediction 
for a random walk in a linear potential (RW-LP model [6]). 

distributions. As Fig. 3 b) demonstrates, the agreement 
is again very satisfactory. 



4. Outlook: Polymer-Membrane Systems 

The success of the phenomenological description 
(1) and (2) motivates its appUcation to more complex 
situations. As an example, we shall now discuss the 
effect of inserting a hydrophilic polymer in our sys- 
tem. Two possible scenarios are illustrated in Fig. 4 
(left): Either the polymer stays squeezed between two 
membranes and deforms them locally, or it spans be- 
tween layers, creating one or several holes. 

We have performed simulations of a membrane 
stack containing a polymer with 100 iJ-beads [7]. In 
free iJ-solution, the polymer exhibits self-avoiding 
walk statistics with the Kuhn length a = 0.4cr. A 
configuration in a membrane stack is shown in Fig. 4 
(right). The polymer induces a hole in one of the ad- 
jacent membranes and collapses into a globule which 
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Fig. 4. Left: Possible ways of incorporating a polymer in a 
membrane stack: A - Squeezed between two membranes; B - 
Spanning several layers and creating a hole. Right: Configuration 
snapshot of a polymer of length A'^ = 100 in the membrane stack. 
The surfaces show the mid-planes of two membranes close to the 
polymer. 

spreads between two layers. 

The situation can be analyzed using scaling argu- 
ments. The entropy penalty for confining a polymer 
of length A'' into a disk of width R± and diameter 

is [8,7] 

Ec - Nia/R^f^'' + N^a^'/iR^Rp i?,, > R^ (3) 
Ec ~ iV3/4(aVi?x-R|)'/^ R\\ < R\i (4) 

with i?| ~ Afi/2^y*^a5/6 £q 3 corresponds to a 
2D compressed (pancake) regime, where the main en- 
tropic penalty comes from the confinement in the R± 
direction (first term), and the polymer behaves like a 
slightly compressed 2D self avoiding walk in the R\\ 
direction (second term). Eq. 4 describes a 3D com- 
pressed regime. 

In the scenario A of Fig. 4, Ep has to be balanced 
with the elastic energy E'gj of fiie membrane stack. 
FromEq. 1, we can derive [7] i?gj ^ 'kKc{R\\ Sy/^^, 
where S is the deformation of the membranes closest 
to the polymer. In the scenario B, the elastic energy 
can be neglected, and the polymer compression energy 
has to be balanced with the energy of pore formation, 
-Epore ~ 27rAi?||. 

Details of the calculation will be given else- 
where [7]. Here we merely sketch some of the results: 
Short chains remain squeezed between membranes 
(scenario A) and assume a 2D compressed (pancake) 
structure. In this regime, the total energy scales hn- 
early with N. Long chains collapse into a 3D com- 
pressed globule, which creates one pore and spans 
exactly two layers (scenario B). Here, the total energy 
scales with N^/^^. In stacks with low compressibil- 



ity (large BS^), 2D compressed states spanning two 
or several layers may appear at intermediate chain 
lengths. 

This quaUtatively explains the observation in the 
simulation: The polymer is in the long chain regime, 
collapses into the 3D compressed state and spreads 
between two layers. The agreement is not yet quantita- 
tive. According to the scaling theory, the transition to 
the long-chain regime is expected at the chain length 
~ A^a-^/^J^^/^, which corresponds to ~ 2000 
with our parameter values. Our chain is much shorter 
The discrepancy can presumably be explained by the 
unknown prefactor in the scaling law for Nc- 

To conclude, the phenomenological models (1) and 
(2) are useful to analyze not only pure membrane 
stacks, but also polymer-membrane systems. It is 
remarkable that both the scaling theory and the simu- 
lations predict that long (non-adsorbing) hydrophilic 
polymers enforce pore formation in the membranes. 
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